//----------------------------------*-C++-*----------------------------------//
/*!
 * \file   Moments_Field.hh
 * \author Jeremy Roberts
 * \date   Jul 6, 2011
 * \brief  
 * \note   Copyright (C) 2011 Jeremy Roberts. 
 */
//---------------------------------------------------------------------------//
// $Rev::                                               $:Rev of last commit
// $Author::                                            $:Author of last commit
// $Date::                                              $:Date of last commit
//---------------------------------------------------------------------------//

#ifndef MOMENTS_FIELD_HH_
#define MOMENTS_FIELD_HH_

#include <cmath>

#include "tran/Cell_Field.hh"
#include "angle/AngleTraits.hh"
#include "../utilities/SP.hh"

// equations
#include "tran/DD_1D_Equations.hh"
#include "tran/DD_2D_Equations.hh"
#include "tran/DD_3D_Equations.hh"

namespace slabtran
{

//===========================================================================//
/*!
 * \class Moments_Field
 * \brief A container of moments and related accessors.
 *
 * Note, the accessors always come in consts/mutable pairs.  Following the
 * approach by Scott Meyers, we implement the const version and have the
 * mutable version call the const version via type casting.  This reduces
 * future maintenance issues.
 *
 * \refs
 * - Scott Meyers, <em>Effective C++</em>, Addision Wesley, 2005.  See in
 *   particular "Item 3".
 *
 * \tparam <Equations> equation type
 * \tparam <D>         dimension
 */
 /*!
 * \example test/testMoments_Field.cc
 *
 * Test of Moments_Field.
 *
 */
//===========================================================================//

template <class Equations, class D>
class Moments_Field : public Cell_Field<double>
{
public:

  /// \name Typedefs
  //\{
  typedef typename util::SP<Moments_Field<Equations, D> >   SP_moments_field;
  typedef Cell_Field<double>                                Base;
  typedef typename Base::value_type                         value_type;
  typedef typename Base::reference                          reference;
  typedef typename Base::const_reference                    const_reference;
  typedef typename Base::pointer                            pointer;
  typedef typename Base::const_pointer                      const_pointer;
  typedef typename Base::iterator                           iterator;
  typedef typename Base::const_iterator                     const_iterator;
  typedef typename Base::size_type                          size_type;
  typedef typename Base::difference_type                    difference_type;
  // convenience type
  typedef typename Base::SP_mesh                            SP_mesh;
  //\}

  /*!
   * \brief Constructor.
   *
   * \param mesh              Mesh we live on.
   * \param legendre_order    Number of moments for the <e>flux expansion</e>.
   * \param initial_value     Initial value for all members of field.
   */
  Moments_Field(SP_mesh mesh, int legendre_order, double initial_value = 0.0);

  /// \name Accessors
  //\{

  /*!
   * \brief Const reference to moment.
   *
   * This is a general interface for all dimensions based on cardinal
   * cell and moment indices.  This is the preferred interface for
   * inner products representing loops over moments.
   *
   * \param     cell    Cardinal spatial index.
   * \param     u       Cell unknown (depends on \em Equations).
   * \param     moment  Cardinal moment index.
   * \return            Const reference to moment value.
   */
  const value_type& operator()(int cell, int u, int moment) const;

  /*!
   * \brief Mutable reference to moment.
   *
   * This is a general interface for all dimensions based on cardinal
   * cell and moment indices.  This mutable reference
   *
   * \param     cell    Cardinal spatial index.
   * \param     u       Cell unknown (depends on \em Equations).
   * \param     moment  Cardinal moment index.
   * \return            Mutable reference to moment value.
   */
  value_type& operator()(int cell, int u, int moment);

  /*!
   * \brief Const reference to moment.
   *
   * This is a general interface for all dimensions.  While the "moment"
   * interface methods below offer an easy way to get moments within a
   * space and angle loop, this operator() makes it easy to use the same
   * operator for all dimensions.  The client is obligated to use proper
   * l and m values.
   *
   * \param s           Cardinal spatial index.
   * \param u           Cell unknown (depends on \em Equations).
   * \param l           Legendre degree.
   * \param m           Legendre order.
   * \return            Const reference to moment value.
   */
  const value_type& operator()(int cell, int u, int l, int m) const;

  /*!
   * \brief Mutable reference to moment.
   *
   * This is a general interface for all dimensions.  While the "moment"
   * interface methods below offer an easy way to get moments within a
   * space and angle loop, this operator() makes it easy to use the same
   * operator for all dimensions.  The client is obligated to use proper
   * l and m values.
   *
   * \param s           cardinal spatial index
   * \param u           cell unknown (depends on \em Equations)
   * \param l           Legendre degree
   * \param m           Legendre order
   * \return            Mutable reference to moment value.
   */
  value_type& operator()(int cell, int u, int l, int m);

  /*!
   * \brief Const reference to moment via ijk index.
   *
   * This is the general interface for 3d.
   *
   * \param     i   X index.
   * \param     j   Y index.
   * \param     k   Z index.
   * \param     u   Cell unknown.
   * \param     l   Spherical harmonic degree.
   * \param     m   Spherical harmonic order.
   * \return        Const reference to moment value.
   */
  const value_type& moment(int i, int j, int k, int u, int l, int m) const;

  /*!
   * \brief Mutable reference to moment via ijk index.
   *
   * This is the general interface for 3d.
   *
   * \param     i   X index.
   * \param     j   Y index.
   * \param     k   Z index.
   * \param     u   Cell unknown.
   * \param     l   Spherical harmonic degree.
   * \param     m   Spherical harmonic order.
   * \return        Mutable reference to moment value.
   */
  value_type& moment(int i, int j, int k, int u, int l, int m);

  /*!
   * \brief Const reference to moment via ijk index.
   *
   * This is the interface for use in 2d computations.  It
   * requires that legal values of \f$(l,m)\f$ for 2d be used; else, the
   * caller is confused!
   *
   * \param     i   X index.
   * \param     j   Y index.
   * \param     k   Z index.
   * \param     u   Cell unknown.
   * \param     l   Spherical harmonic degree.
   * \param     m   Spherical harmonic order.
   * \return        Const reference to moment value.
   */
  const value_type& moment(int i, int j, int u, int l, int m) const;

  /*!
   * \brief Mutable reference to moment via ijk index.
   *
   * This is the interface for use in 2d computations.  It
   * requires that legal values of \f$(l,m)\f$ for 2d be used; else, the
   * caller is confused!
   *
   * \param     i   X index.
   * \param     j   Y index.
   * \param     k   Z index.
   * \param     u   Cell unknown.
   * \param     l   Spherical harmonic degree.
   * \param     m   Spherical harmonic order.
   * \return        Mutable reference to moment value.
   */
  value_type& moment(int i, int j, int u, int l, int m);

  /*!
  * \brief Const reference to moment via ijk index.
  *
  * This is the interface for use in 1d computations.
  *
  * \param     i   X index.
  * \param     u   Cell unknown.
  * \param     l   Legendre order.
  * \return        Const reference to moment value.
  */
  const value_type& moment(int i, int u, int l) const;


   /*!
   * \brief Mutable reference to moment via ijk index.
   *
   * This is the interface for use in 1d computations.
   *
   * \param     i   X index.
   * \param     u   Cell unknown.
   * \param     l   Legendre order.
   * \return        Mutable reference to moment value.
   */
  value_type& moment(int i, int u, int l);

  /*!
   *  \brief Return Legendre order
   */
  int legendre_order() const { return d_legendre_order; };

  /*!
   *  \brief Return number of spherical harmonic moments per unknown.
   */
  int number_moments() const { return d_number_moments; };

  //\}

private:

  /// \name Data
  //\{

  /// Legendre order (as in \f$N\f$ of \f$P_N\f$).
  const int d_legendre_order;

  /// Number of spherical harmonic flux moments per unknown.
  const int d_number_moments;

  //\}

  /// \name Implementation
  //\{

  /// Set the m indices based on dimension
  void set_m_index();

  //\}

};


} // end namespace slabtran

//---------------------------------------------------------------------------//
// INLINE FUNCTIONS
//---------------------------------------------------------------------------//

#include "Moments_Field.i.hh"

#endif /* MOMENTS_FIELD_HH_ */

//---------------------------------------------------------------------------//
//              end of Moments_Field.hh
//---------------------------------------------------------------------------//
